ggplot(aes(x=gasmethane))+
geom_density()+
geom_area(aes(x=stage(gasmethane,after_stat=oob_censor(x,c(69.9,94.9)))),
stat="density",fill=brewer.pal(8,"Dark2")[1],alpha=.5)+
labs(x="Methane as percent of natural gas composition",y="Density of responses")+
theme_classic()+
theme(axis.text.y=element_blank(),
axis.ticks.y=element_blank())
## number of responses in the "correct" range
nrow(d[d$gasmethane>69.9&d$gasmethane<95,])/nrow(d)
ggsave(paste0(figures,"descriptive_methanegascomp.pdf"),width=6.5,height=6.5,
ggarrange(fig1+theme(legend.position="left"),gascompfig,nrow=2,ncol=1))
## Fig. 2: knowledge, perceptions, concern about methane ####
d<-d%>%
mutate(gascomposition.binary=ifelse(gascomposition=="Methane",1,0))
### responses to change in methane in atmosphere question ####
methanechangetab<-data.frame(table(d$methanechange,d$country))
methanechangen<-methanechangetab%>%
group_by(Var2)%>%
summarise(n=sum(Freq))
methanechangetab<-methanechangetab%>%left_join(methanechangen,by="Var2")%>%
mutate(prop=Freq*100/n)%>%
select(Var2,Var1,prop)%>%
mutate(Var1=factor(Var1,levels=c("stayed same","decreased","increased")))
fig2<-methanechangetab%>%
ggplot(aes(x=prop,y=Var2,fill=Var1))+
theme_bw()+
geom_bar(stat="identity",position="stack")+
# scale_fill_brewer(palette="Dark2",guide=guide_legend(reverse=TRUE))+
scale_fill_manual(values=brewer.pal(8,"Dark2")[2:4],guide=guide_legend(reverse=TRUE))+
labs(x="Percent of respondents",y="",fill="Has methane in the\natmosphere increased,\ndecreased, or stayed\nthe same?")+
theme(legend.position="left")
fig2
### question about methane as a problem ####
probtab<-data.frame(table(d$country,d$methaneproblem))
probn<-probtab%>%
group_by(Var1)%>%
summarise(n=sum(Freq))%>%
ungroup()
probtab<-left_join(probtab,probn,by="Var1")%>%
mutate(prop=Freq*100/n)%>%
select(Var1,Var2,prop)
probfig<-probtab%>%
ggplot(aes(x=prop,y=Var1,fill=Var2))+
theme_bw()+
geom_bar(stat="identity",position="stack")+
# scale_fill_brewer(palette="Dark2",guide=guide_legend(reverse=TRUE))+
scale_fill_manual(values=brewer.pal(8,"Dark2")[c(8,2:4)],guide=guide_legend(reverse=TRUE))+
labs(x="Percent of respondents",y="",fill="How much of a problem\nis methane\nfor the climate?")+
theme(legend.position="left")
probfig
# save figure as two panels
ggsave(paste0(figures,"descriptive_methaneknowledge_beliefs.pdf"),width=6.5,height=6.5,
ggarrange(fig2+labs(x=""),probfig,
ncol=1,nrow=2))
# Fig 3: policy questions that precede our experiment ####
## collapse methane fee questions (cost info doesn't make much of a difference, so note in paper that this marginalizes across random assignment of cost and benefit info)
fee<-d%>%
select(importregulations_2a,importregulations_3a,importregulations_4a)%>%
pivot_longer(everything(),names_to="var",values_to="val")%>%
filter(!is.na(val))
feesum<-data.frame(prop.table(table(fee$val))*100)%>%mutate(var="Impose fee on\nimported oil and gas")
importregs<-data.frame(prop.table(table(d$importregulations_1a))*100)%>%mutate(var="Regulate imported gas")
emissionregs<-data.frame(prop.table(table(d$methaneemissiona))*100)%>%mutate(var="Regulate methane leaks")
regcompanies<-data.frame(prop.table(table(d$financialincentive))*100)%>%mutate(var="Regulate companies")
techupdate<-data.frame(prop.table(table(d$techupdatea))*100)%>%mutate(var="Require equipment upgrades")
regfig<-do.call(rbind,list(feesum,importregs,emissionregs,techupdate))%>%
mutate(response=case_when(Var1==0~"Strongly oppose",
Var1==1~"Oppose",
Var1==2~"Support",
Var1==3~"Strongly support"),
response=factor(response,levels=c("Strongly support","Support","Oppose","Strongly oppose")))%>%
ggplot(aes(x=Freq,y=var,fill=response))+
theme_bw()+
geom_bar(stat="identity",position="stack")+
scale_fill_brewer(palette="Dark2",guide=guide_legend(reverse=TRUE))+
labs(x="Percent of respondents",y="Policy proposal",fill="Response")
regfig
ggsave(file=paste0(figures,"descriptive_policysupport.pdf"),width=6.5,height=3.5,regfig)
# Experimental analysis  ####
## cost sensitivity regressions ####
# set up IPW weights
sum(table(d$z.cost))
w <- weightit(z.cost~block,data=d,estimand="ATE",method="ps")
d$w.cost<-w$weights
d$z.cost<-relevel(factor(d$z.cost),ref="no cost info")
h1mod<-lm_robust(dv.methanelaw.binary~z.cost+country+ideo+gender+age+education_recode+income_recode,weights=w.cost,data=d)
## apply Benjamini-Hochberg correction
bh.h1mod <- round(p.adjust(h1mod$p.value[grepl("z",names(h1mod$p.value))],"BH"),2)
## analysis with ordinal outcome as continuous variable
h1mod.cont<-lm_robust(dv.methanelaw~z.cost+country+ideo+gender+age+education_recode+income_recode,weights=w.cost,data=d)
bh.h1mod.cont <- round(p.adjust(h1mod.cont$p.value[grepl("z",names(h1mod$p.value))],"BH"),2)
## Table S6: cost results
texreg(list(h1mod,h1mod.cont),include.ci = FALSE,custom.model.names=c("LPM Model","Linear model with ordinal DV"),
caption=paste("Effect of cost information on policy support. Benjamini-Hochberg adjusted p-values for cost coefficients:", bh.h1mod[1],  bh.h1mod[2], "(LPM)",bh.h1mod.cont[1],bh.h1mod.cont[2],"(Continuous DV)", sep = " "),
file=paste0(tables,"costresults.tex"),
label="tab:costresults",float.pos="!htpb")
## benefit frames regressions ####
w<-weightit(z.frame~block,estimand="ATE",method="ps",data=d)
d$w.frame<-w$weights
d$z.frame<-relevel(factor(d$z.frame),ref="no frame")
h2mod<-lm_robust(dv.methanelaw.binary~z.frame+country+ideo+gender+age+education_recode+income_recode,weights=w.frame,data=d)
## extract Benjamini-Hochberg corrected p values
bh.h2mod <- round(p.adjust(h2mod$p.value[grepl("z.frame",names(h2mod$p.value))],"BH"),2)
## continuous DV
h2mod.cont<-lm_robust(dv.methanelaw~z.frame+country+ideo+gender+age+education_recode+income_recode,weights=w.frame,data=d)
## extract Benjamini-Hochberg corrected p values
bh.h2mod.cont <- round(p.adjust(h2mod.cont$p.value[grepl("z.frame",names(h2mod$p.value))],"BH"),2)
## table S7: framing results
texreg(list(h2mod,h2mod.cont),include.ci=FALSE,custom.model.names = c("LPM Model","Linear model with ordinal DV"),
caption=paste("Effects of benefits framing on policy support. Benjamini-Hochberg adjusted p-values for benefit coefficients:", bh.h2mod[1], bh.h2mod[2], bh.h2mod[3], "(LPM model)",
bh.h2mod.cont[1], bh.h2mod.cont[2], bh.h2mod.cont[3], "(Linear model with ordinal DV)",sep = " "),
file=paste0(tables,"frameresults.tex"),
label="tab:frameresults",float.pos="!htpb")
## cost x benefit regressions ####
## code interacted treatment groups
table(d$z.cost)
table(d$z.frame)
d<-d%>%
mutate(z=case_when(z.cost=="no cost info"&z.frame=="no frame"~"control",
z.cost=="high"&z.frame=="climate"~"increase.climate",
z.cost=="low"&z.frame=="climate"~"decrease.climate",
z.cost=="high"&z.frame=="energysecurity"~"increase.security",
z.cost=="low"&z.frame=="energysecurity"~"decrease.security",
z.cost=="high"&z.frame=="health"~"increase.health",
z.cost=="low"&z.frame=="health"~"decrease.health"))
w<- weightit(z~block,data=d,estimand="ATE",method="ps")
d$w.full<-w$weights
h3amod<-lm_robust(dv.methanelaw.binary~z+country+ideo+gender+age+education_recode+income_recode,
data=d,weights=w.full)
h3amod.cont<-lm_robust(dv.methanelaw~z+country+ideo+gender+age+education_recode+income_recode,
data=d,weights=w.full)
## table S8: cost sensitivity of benefits framing
texreg(list(h3amod,h3amod.cont),
custom.model.names=c("LPM Model","Linear model with ordinal DV"),
custom.coef.names = c("(Intercept)",
"Cost.decrease:Benefit.climate",
"Cost.decrease:Benefit.health",
"Cost.decrease:Benefit.security",
"Cost.increase:Benefit.climate",
"Cost.increase:Benefit.health",
"Cost.increase:Benefit.security",
"Germany","Italy","Poland","Left","Right","Male",
"Age","Secondary degree or lower","Income"),
include.ci=FALSE,caption="Cost sensitivity of benefit framing",
file=paste0(tables,"results_interactions.tex"),
label="tab:interactions",float.pos="!htpb")
## wald test for cost-sensitivity of response to benefit frames ####
waldclimate<-as.data.frame(linearHypothesis(h3amod,"zdecrease.climate=zincrease.climate"))[2,]%>%
mutate(`Benefit frame`="Climate")%>%
select(`Benefit frame`,`Chisq`,`Pr(>Chisq)`)
waldhealth<-as.data.frame(linearHypothesis(h3amod,"zdecrease.health=zincrease.health"))[2,]%>%
mutate(`Benefit frame`="Health")%>%
select(`Benefit frame`,`Chisq`,`Pr(>Chisq)`) ## this one is significant, but in the wrong direction. (ie, effect in the cost-increasing condition is positive, but null. effect in cost-decreasing condition is negative, and significant)
waldsecurity<-as.data.frame(linearHypothesis(h3amod,"zdecrease.security=zincrease.security"))[2,]%>%
mutate(`Benefit frame`="Security")%>%
select(`Benefit frame`,`Chisq`,`Pr(>Chisq)`)
wald.tests<-do.call(rbind,list(waldclimate,waldhealth,waldsecurity))
## table S9: interactions between cost and benefit frames ####
print.xtable(xtable(wald.tests,caption="Wald tests for interaction effects between cost and benefit frames",
label="tab:waldtests"),include.rownames = FALSE,file=paste0(tables,"results_waldtests.tex"))
linearHypothesis(h3amod.cont,"zdecrease.climate=zincrease.climate")
linearHypothesis(h3amod.cont,"zdecrease.health=zincrease.health")
linearHypothesis(h3amod.cont,"zdecrease.security=zincrease.security")
## 2-panel plot with marginalized and fully interacted effects ####
frameresults<-tidy(h2mod)
costresults<-tidy(h1mod)
intresults<-tidy(h3amod)
frameresults.cont<-tidy(h2mod.cont)
costresults.cont<-tidy(h1mod.cont)
intresults.cont<-tidy(h3amod.cont)
results<-rbind(frameresults,costresults)
results.cont<-rbind(frameresults.cont,costresults.cont)
mainplot<-results%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.64*std.error,
ci.low=estimate-1.64*std.error,
ci.high95=estimate+1.96*std.error,
ci.low95=estimate-1.96*std.error)%>%
ggplot(aes(x=estimate,y=term))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high))+
geom_linerange(aes(xmin=ci.low95,xmax=ci.high95),lwd=.15)+
geom_vline(xintercept=0,lty="dotted")+
theme_classic()+
scale_y_discrete(
labels = c("z.costlow" = "Lower costs",
"z.costhigh" = "Higher costs",
"z.framehealth" = "Health",
"z.frameclimate" = "Climate change",
"z.frameenergysecurity" = "Energy security"))+
labs(x="",y="")+
coord_cartesian(xlim = c(-0.1,0.05)) +
theme(axis.text.x = element_text(size=11),
axis.text.y = element_text(size=11))
mainplot
mainplot.cont<-results.cont%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.64*std.error,
ci.low=estimate-1.64*std.error)%>%
ggplot(aes(x=estimate,y=term))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high))+
geom_vline(xintercept=0,lty="dotted")+
theme_classic()+
scale_y_discrete(
labels = c("z.costlow" = "Lower costs",
"z.costhigh" = "Higher costs",
"z.framehealth" = "Health",
"z.frameclimate" = "Climate change",
"z.frameenergysecurity" = "Energy security"))+
labs(x="",y="")+
coord_cartesian(xlim = c(-0.2,0.1)) +
theme(axis.text.x = element_text(size=11),
axis.text.y = element_text(size=11))
mainplot.cont
intplot<-intresults%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.64*std.error,
ci.low=estimate-1.64*std.error,
ci.high95=estimate+1.96*std.error,
ci.low95=estimate-1.96*std.error)%>%
ggplot(aes(x=estimate,y=term))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high))+
geom_linerange(aes(xmin=ci.low95,xmax=ci.high95),lwd=.15)+
geom_vline(xintercept=0,lty="dotted")+
scale_y_discrete(limits=c("zincrease.climate","zdecrease.climate","zincrease.security","zdecrease.security","zincrease.health","zdecrease.health"),
labels = c("zincrease.climate"= "Higher cost,\nclimate change","zdecrease.climate"="Lower cost,\nclimate change",
"zincrease.security"="Higher cost,\nenergy security","zdecrease.security"="Lower cost,\nenergy security",
"zincrease.health"="Higher cost,\nhealth","zdecrease.health"="Lower cost,\nhealth")
) +
labs(x="",y="")+
coord_cartesian(xlim = c(-0.1,0.05)) +
theme_classic()+
theme(axis.text.x = element_text(size=11),
axis.text.y = element_text(size=11))
intplot
intplot.cont<-intresults.cont%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.64*std.error,
ci.low=estimate-1.64*std.error)%>%
ggplot(aes(x=estimate,y=term))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high))+
geom_vline(xintercept=0,lty="dotted")+
scale_y_discrete(limits=c("zincrease.climate","zdecrease.climate","zincrease.security","zdecrease.security","zincrease.health","zdecrease.health"),
labels = c("zincrease.climate"= "Higher cost,\nclimate change","zdecrease.climate"="Lower cost,\nclimate change",
"zincrease.security"="Higher cost,\nenergy security","zdecrease.security"="Lower cost,\nenergy security",
"zincrease.health"="Higher cost,\nhealth","zdecrease.health"="Lower cost,\nhealth")
) +
labs(x="",y="")+
coord_cartesian(xlim = c(-0.2,0.1)) +
theme_classic()+
theme(axis.text.x = element_text(size=11),
axis.text.y = element_text(size=11))
intplot.cont
## figure 4: combined results plot ####
results.combined<-annotate_figure(ggarrange(intplot,mainplot),
left=text_grob("Treatment",size=11,rot=90),
bottom = text_grob("Effect on policy support, relative to control", size = 11))
ggsave(file=paste0(figures,"results_combined.pdf"),height=3.5,width=6.5,results.combined)
## figure S2: combined results plot with ordinal DV
results.combined.cont<-annotate_figure(ggarrange(intplot.cont,mainplot.cont),
left=text_grob("Treatment",size=11,rot=90),
bottom = text_grob("Effect on policy support, relative to control", size = 11))
ggsave(file=paste0(figures,"results_combined_continuousDV.pdf"),height=3.5,width=6.5,
results.combined.cont)
# minimum detectable effect ####
## standard error of the average treatment effect for the low-cost treatment:
## pg. 61 from Gerber and Green 2012 (eqn. 3.6)
## var Y0:
d2<-d%>%
filter(z.cost=="no cost info"|z.cost=="high")%>%
mutate(d.diff=case_when(z.cost=="no cost info"~dv.methanelaw.binary-mean(d$dv.methanelaw.binary[d$z.cost=="no cost info"]),
z.cost=="high"~dv.methanelaw.binary-mean(d$dv.methanelaw.binary[d$z.cost=="high"])))
se_y0_hat<-sqrt(sum((d2$d.diff[d2$z.cost=="no cost info"])^2)/(nrow(d2[d2$z.cost=="no cost info",])-1))
se_y1_hat<-sqrt(sum((d2$d.diff[d2$z.cost=="high"])^2)/(nrow(d2[d2$z.cost=="high",])-1))
# sd_y_0_hat <- sd(d$dv.methanelaw.binary[d$z.cost == "no cost info"])
# sd_y_1_hat <- sd(d$dv.methanelaw.binary[d$z.cost == "high"])
# sd_y_2_hat <- sd(d$dv.methanelaw.binary[d$z.cost == "low"])
se_hat_full_study <- sqrt( se_y0_hat^2 / nrow(d2[d2$z.cost=="no cost info",])  + se_y1_hat^2 / nrow(d2[d2$z.cost=="high",]))
mde_hat_full_study <- 2.8 * se_hat_full_study
# results by knowledge ####
table(d$gascomposition.binary)
table(d$z.cost)
## cost regressions by knowledge ####
d$knowledge<-factor(d$gascomposition.binary,levels=c(0,1),labels=c("low","high"))
w <- weightit(z.cost~block,data=d[d$gascomposition.binary==1,],estimand="ATE",method="ps")
d[d$gascomposition.binary==1,"w.cost.byknowledge"]<-w$weights
w<-weightit(z.cost~block,data=d[d$gascomposition.binary==0,],estimand="ATE",method="ps")
d[d$gascomposition.binary==0,"w.cost.byknowledge"]<-w$weights
h1mod.know1<-tidy(lm_robust(dv.methanelaw.binary~z.cost+country+ideo+gender+age+education_recode+income_recode,
weights=w.cost.byknowledge,data=subset(d,gascomposition.binary==1)))%>%
mutate(knowledge="high")
h1mod.know0<-tidy(lm_robust(dv.methanelaw.binary~z.cost+country+ideo+gender+age+education_recode+income_recode,
weights=w.cost.byknowledge,data=subset(d,gascomposition.binary==0)))%>%
mutate(knowledge="low")
costint<-tidy(lm_robust(dv.methanelaw.binary~z.cost*knowledge+country+ideo+gender+age+education_recode+income_recode,
weights=w.cost,data=d))%>%
filter(grepl(":",term))%>%
mutate(sig=ifelse(p.value<0.05,1,0))%>%
select(term,p.value,sig)## neither effect is sig
h1mod.knowledge<-rbind(h1mod.know1,h1mod.know0)
## benefit frames regressions by knowledge ####
w<-weightit(z.frame~block,estimand="ATE",method="ps",data=d[d$gascomposition.binary==1,])
d[d$gascomposition.binary==1,"w.frame.know"]<-w$weights
w<-weightit(z.frame~block,estimand="ATE",method="ps",data=d[d$gascomposition.binary==0,])
d[d$gascomposition.binary==0,"w.frame.know"]<-w$weights
frameresults.know1<-tidy(lm_robust(dv.methanelaw.binary~z.frame+country+ideo+gender+age+education_recode+income_recode,
weights=w.frame.know,data=subset(d,gascomposition.binary==1)))%>%
mutate(knowledge="high")
frameresults.know0<-tidy(lm_robust(dv.methanelaw.binary~z.frame+country+ideo+gender+age+education_recode+income_recode,
weights=w.frame.know,data=subset(d,gascomposition.binary==0)))%>%
mutate(knowledge="low")
h2mod.knowledge<-rbind(frameresults.know1,frameresults.know0)
frameint<-tidy(lm_robust(dv.methanelaw.binary~z.frame*knowledge+country+ideo+gender+age+education_recode+income_recode,
weights=w.frame,data=d))%>%
mutate(sig=ifelse(p.value<0.05,1,0))%>%
filter(grepl(":",term))%>%
select(term,p.value,sig) ## all null
## interaction effects by knowledge ####
w<- weightit(z~block,data=d[d$gascomposition.binary==1,],estimand="ATE",method="ps")
d[d$gascomposition.binary==1,"w.full.know"]<-w$weights
w<- weightit(z~block,data=d[d$gascomposition.binary==0,],estimand="ATE",method="ps")
d[d$gascomposition.binary==0,"w.full.know"]<-w$weights
h3amod.know1<-tidy(lm_robust(dv.methanelaw.binary~z+country+ideo+gender+age+education_recode+income_recode,
data=subset(d,gascomposition.binary==1),weights=w.full.know))%>%
mutate(knowledge="high")
h3amod.know0<-tidy(lm_robust(dv.methanelaw.binary~z+country+ideo+gender+age+education_recode+income_recode,
data=subset(d,gascomposition.binary==0),weights=w.full.know))%>%
mutate(knowledge="low")
h3amod.know<-rbind(h3amod.know0,h3amod.know1)
fullint<-tidy(lm_robust(dv.methanelaw.binary~z*knowledge+country+ideo+gender+age+education_recode+income_recode,
weights=w.frame,data=d))%>%
filter(grepl(":",term))%>%
mutate(sig=ifelse(p.value<0.05,1,0))%>%
select(term,p.value,sig) ## increase.health is only significant difference
## table S10 of results by knowledge ####
resultstab.know<-do.call(rbind,list(h3amod.know,h1mod.knowledge,h2mod.knowledge))%>%
filter(grepl("z",term))
resultstab.know$`BH p value` <- round(p.adjust(resultstab.know$p.value,"BH"),2)
resultstab.know<-resultstab.know%>%
select(term,estimate,std.error,statistic,p.value,`BH p value`,knowledge)%>%
mutate(term=case_when(term=="zincrease.security"~"Higher cost, energy security",
term=="zdecrease.security"~"Lower cost, energy security",
term=="zincrease.climate"~"Higher cost, climate change",
term=="zdecrease.climate"~"Lower cost, climate change",
term=="zincrease.health"~"Higher cost, health",
term=="zdecrease.health"~"Lower cost, health",
term=="z.costlow"~ "Lower costs",
term=="z.costhigh"~"Higher costs",
term=="z.framehealth"~"Health",
term=="z.frameclimate"~"Climate change",
term=="z.frameenergysecurity"~"Energy security"))
print.xtable(xtable(resultstab.know,caption="{\\bf Experimental effects by knowledge:} The table shows the results from OLS regressions of our treatment indicators for the cost, benefit, and cost x benefit treatments on split samples of respondents defined by their knowledge about methane in the atmosphere. The penultimate column shows p-values after applying the Benjamini-Hochberg correction for multiple hypothesis testing.",
label="tab:knowledgeresults"),include.rownames = FALSE,file=paste0(tables,"results_knowledgeresults.tex"))
## Fig. S3: 2-panel effects plot by knowledge ####
knowplot.right<-rbind(h1mod.knowledge,h2mod.knowledge)%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.96*std.error,
ci.low=estimate-1.96*std.error)%>%
ggplot(aes(x=estimate,y=term,shape=knowledge))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high),position=position_dodge(width=.2))+
geom_vline(xintercept=0,lty="dotted")+
theme_classic()+
scale_y_discrete(
labels = c("z.costlow" = "Lower costs",
"z.costhigh" = "Higher costs",
"z.framehealth" = "Health",
"z.frameclimate" = "Climate change",
"z.frameenergysecurity" = "Energy security"))+
labs(x="",y="")+
scale_x_continuous(limits=c(-.15,.15))+
# coord_cartesian(xlim = c(-0.1,0.1)) +
theme(axis.text.x = element_text(size=11),
axis.text.y = element_text(size=11))
knowplot.right
knowplot.left<-h3amod.know%>%
filter(grepl("z.",term))%>%
mutate(ci.high=estimate+1.96*std.error,
ci.low=estimate-1.96*std.error)%>%
ggplot(aes(x=estimate,y=term,shape=knowledge))+
geom_pointrange(aes(xmin=ci.low,xmax=ci.high),position=position_dodge(width=.2))+
geom_vline(xintercept=0,lty="dotted")+
scale_y_discrete(
limits = c("zincrease.security","zdecrease.security",
"zincrease.climate","zdecrease.climate",
"zincrease.health","zdecrease.health"),
labels = c("Higher cost,\nenergy security","Lower cost,\nenergy security",
"Higher cost,\nclimate change","Lower cost,\nclimate change",
"Higher cost,\nhealth","Lower cost,\nhealth")
) +
labs(x="",y="")+
scale_x_continuous(limits=c(-.15,.15))+
# coord_cartesian(xlim = c(-0.1,0.1)) +
theme(legend.position = "bottom",
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.background = element_blank(),
axis.line.x = element_line(colour = 'black', size=0.5, linetype='solid'),
axis.text.x = element_text(size=11),
axis.title.x = element_text(size=13),
axis.line.y = element_line(colour = 'black', size=0.5, linetype='solid'),
axis.text.y = element_text(size=11),
axis.title.y = element_text(size=13))
knowplot.left
results.combined.byknowledge<-annotate_figure(ggarrange(knowplot.left,knowplot.right,common.legend=TRUE,
legend="bottom"),
left=text_grob("Treatment",size=11,rot=90),
bottom = text_grob("Effect on policy support, relative to control", size = 11))
ggsave(file=paste0(figures,"results_combined_byknowledge.pdf"),height=3.5,width=6.5,results.combined.byknowledge)
# exploratory hypotheses ####
## cost x covariate interactions ####
eh1amod<-lm_robust(dv.methanelaw.binary~z.cost*country+ideo+gender+age+education_recode+income_recode,
weights=w.cost,data=d)
eh1bmod<-lm_robust(dv.methanelaw.binary~z.cost*ideo+country+gender+age+education_recode+income_recode,
weights=w.cost,data=d)
eh1cmod<-lm_robust(dv.methanelaw.binary~z.cost*age+country+ideo+gender+education_recode+income_recode,
weights=w.cost,data=d)
eh1dmod<-lm_robust(dv.methanelaw.binary~z.cost*education_recode+ideo+gender+age,#+income+education,
weights=w.cost,data=d)
eh1emod<-lm_robust(dv.methanelaw.binary~z.cost*income_recode+ideo+gender+age,#+income+education,
weights=w.cost,data=d)
## table S11: conditional effects of cost information ####
texreg(list(eh1amod,eh1bmod, eh1cmod,eh1dmod,eh1emod),
custom.model.names = c("x countries","x ideologies","x ages",
"x education","x incomes"),
include.ci=FALSE,caption="Conditional effects of cost information",
custom.coef.names=c("(Intercept)","High cost","Low cost","Germany","Italy","Poland",
"Left","Right","Male","Age",
"Secondary ed or less","Income",
"High cost x Germany","Low cost x Germany",
"High cost x Italy","Low cost x Italy",
"High cost x Poland","Low cost x Poland",
"High cost x Left","Low cost x Left",
"High cost x Right","Low cost x Right",
"High cost x Age","Low cost x Age",
"High cost x Ed","Low cost x Ed",
"High cost x Income","Low cost x Income"),
fontsize="footnotesize",single.row=TRUE,
file=paste0(tables,"results_costinteractions.tex"),
label="tab:costinteractions")
## benefit x covariate interactions ####
eh2amod<-lm_robust(dv.methanelaw.binary~z.frame*country+ideo+gender+age+education_recode+income_recode,
weights=w.frame,data=d)
eh2bmod<-lm_robust(dv.methanelaw.binary~z.frame*ideo+country+gender+age+education_recode+income_recode,
weights=w.frame,data=d)
eh2cmod<-lm_robust(dv.methanelaw.binary~z.frame*age+country+ideo+gender+education_recode+income_recode,
weights=w.frame,data=d)
eh2dmod<-lm_robust(dv.methanelaw.binary~z.frame*education_recode+country+ideo+gender+age+income_recode,weights=w.frame,data=d) ##  add this and next regression in when we have education and income data coded.
eh2emod<-lm_robust(dv.methanelaw.binary~z.frame*income_recode+country+ideo+gender+education_recode+age,weights=w.frame,data=d)
## table S12: conditional effects of benefit framing ####
texreg(list(eh2amod,eh2bmod,eh2cmod,eh2dmod,eh2emod
),
custom.model.names = c("x countries","x ideologies","x ages","x education","x incomes"
),
custom.coef.names=c("(Intercept)","Climate frame","Security frame","Health frame","Germany","Italy","Poland",
"Left","Right","Male","Age",
"Secondary ed or less","Income",
"Climate frame x Germany","Security frame x Germany",
"Health frame x Germany",
"Climate frame x Italy","Security frame x Italy",
"Health frame x Italy",
"Climate frame x Poland","Security frame x Poland",
"Health frame x Poland",
"Climate frame x Left","Security frame x Left",
"Health frame x Left",
"Climate frame x Right","Security frame x Right",
"Health frame x Left",
"Climate frame x Age","Security frame x Age",
"Health frame x Age",
"Climate frame x Ed","Security frame x Ed",
"Health frame x Ed",
"Climate frame x Income","Security frame x Income",
"Health frame x"),
include.ci=FALSE,caption="Conditional effects of benefit framing",
fontsize="footnotesize",single.row=TRUE,
file=paste0(tables,"results_frameinteractions.tex"),
label="tab:frameinteractions",float.pos="!htpb")
## conditional effects of cost-sensitivity of responses to benefit frames ####
eh3amod<-lm_robust(dv.methanelaw.binary~z*country+ideo+gender+age+education_recode+income_recode,
weights=w.full,data=d)
eh3bmod<-lm_robust(dv.methanelaw.binary~z*ideo+country+gender+age+education_recode+income_recode,
weights=w.full,data=d)
eh3cmod<-lm_robust(dv.methanelaw.binary~z*age+country+ideo+gender+education_recode+income_recode,
weights=w.full,data=d)
eh3dmod<-lm_robust(dv.methanelaw.binary~z*education_recode+country+ideo+gender+age+income_recode,weights=w.full,data=d)
eh3emod<-lm_robust(dv.methanelaw.binary~z*income_recode+country+ideo+gender+education_recode+age,weights=w.full,data=d)
## table S13: conditional effects of cost sensitivity of benefit frames ####
texreg(list(eh3amod,eh3bmod,eh3cmod,eh3dmod,eh3emod
),
custom.model.names = c("x countries","x ideologies","x ages","x education","x incomes"
),caption="Conditional effects of cost sensitivity of benefit frames",
float.pos="!htpb",
include.ci=FALSE,fontsize="tiny",single.row=TRUE,
custom.coef.names=c("Intercept","Low cost x climate","Low cost x health","Low cost x security",
"High cost x climate","High cost x health","High cost x security",
"Germany","Italy","Poland","Left","Right","Male","Age","Secondary ed or less","Income",
"Low cost x climate x Germany","Low cost x health x Germany","Low cost x security x Germany",
"High cost x climate x Germany","High cost x health x Germany","High cost x security x Germany",
"Low cost x climate x Italy","Low cost x health x Italy","Low cost x security x Italy",
"High cost x climate x Italy","High cost x health x Italy","High cost x security x Italy",
"Low cost x climate x Poland","Low cost x health x Poland","Low cost x security x Poland",
"High cost x climate x Poland","High cost x health x Poland","High cost x security x Poland",
"Low cost x climate x Left","Low cost x health x Left","Low cost x security x Left",
"High cost x climate x Left","High cost x health x Left","High cost x security x Left",
"Low cost x climate x Right","Low cost x health x Right","Low cost x security x Right",
"High cost x climate x Right","High cost x health x Right","High cost x security x Right",
"Low cost x climate x age","Low cost x health x age","Low cost x security x age",
"High cost x climate x age","High cost x health x age","High cost x security x age",
"Low cost x climate x ed","Low cost x health x ed","Low cost x security x ed",
"High cost x climate x ed","High cost x health x ed","High cost x security x ed",
"Low cost x climate x income","Low cost x health x income","Low cost x security x income",
"High cost x climate x income","High cost x health x income","High cost x security x income"),
file=paste0(tables,"results_fullinteractions.tex"),label="tab:fullinteractions")
